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A topologically ordered phase on a torus possesses degenerate ground states that transform nontrivially under 
the modular transformations of the torus, generated by Dehn twists. Representation of modular transformations 
on the ground states (modular matrices) characterizes the topological order. We show that the modular matrices 
can be numerically measured as the non-Abelian Berry phase of adiabatic deformations of the lattice model 
placed on a torus. We apply this method to the example of a gauged p^ + ipy superconductor, and show that the 
result is consistent with the topological quantum field theory descriptions. 


When a two-dimensional topologically ordered phase’ in¬ 
habits a manifold, it possesses a robust ground state degener¬ 
acy^ dependent on the topology of the manifold. The topo¬ 
logical ground state degeneracy can not be lifted by any lo¬ 
cal perturbation up to exponentially small corrections. The 
topological order (or the long-range entanglement’) in these 
degenerate ground states can be characterized by the unitary 
transformations among them induced by the diffeomorphisms 
of the manifold, which can be described by topological quan¬ 
tum field theories(TQFT) at low energy. Particularly interest¬ 
ing are the “large” diffeomorphisms that are not continuously 
connected to the identity, whose equivalence classes under 
the identity component form the mapping class group (MCG) 
of the manifold. It is widely believed that the (projective) 
representations of the MCG on the topologically-degenerate 
ground states, defined as the non-Abelian Berry phases of the 
ground states under adiabatic deformations of the manifold'’, 
can be used to characterize the topological phase’’^'. 

In two dimensions, the MCG of a torus is the modular group 
SL(2, Z), generated by the Dehn twist and tt/2 rotation. The 
matrix representations of the two generators of SL(2,Z) on 
the topological ground states, called the modular S and T ma¬ 
trices, are believed to completely characterize the topological 
phase. For example, one can read off the quantum dimensions 
of quasiparticles, fusion rules and braiding statistics from the 
modular matrices. It is therefore an important question to de¬ 
termine the modular matrices from microscopic Hamiltonians 
or ground state wavefunctions in order to identify the topo¬ 
logical order. In continuous space (e.g. quantum Hall sys¬ 
tems) Dehn twists can be realized by adiabatically deform¬ 
ing the metric and was exploited to find the Hall viscosity of 
quantum fluids*"’”. Recently a number of proposals to extract 
the modular matrices numerically have appeared”"’*, mostly 
from entanglement measurement or ground state wavefunc- 
tion overlap. 

While these proposals have been proven successful in nu¬ 
merical applications'-’’^’’”"", the entanglement-based meth¬ 
ods'-’’* mostly work in a cylinder geometry, and the methods 
based on wavefunction overlap”’’” do not actually calculate 
the adiabatic Berry phases. In this work we show that for lat¬ 
tice models on a torus, Dehn twists, or “shear” deformation of 
the system can be directly implemented by adiabatically re¬ 
connecting the bonds of the lattice (see Fig. 1(a) for an illus¬ 
tration of the square lattice). Generally in a finite-size system 
(e.g. a L X L torus ), the resulting adiabatic Berry phase 0 



FIG. 1: (Colored online.) (a) One-step twist by smoothly weaken 
the horizontal red bonds and strengthen the diagonal blue bonds, (b) 
Four flux configurations (boundary conditions) labeled by <E>j,). 
Red bounds are with minus signs. The configurations (0, tt) and 
(tt, 0) are connected by a full Dehn twist along the direction shown 
in (a). 


has the following dependence on system size: 

e = a + + 0{L-^). ( 1 ) 

Here a is a universal quantity, determined by the chiral cen¬ 
tral charge c_ mod 24 and topological twists of quasiparticles. 
The Dehn twist T matrix can be obtained by measuring this 
non-Abelian Berry phase in the ground state manifold. Once 
the Dehn twists along both the meridian (Ty) and longitudinal 
(Tx) cycles of the torus are found, the S matrix can be derived 
from the composition of Dehn twists as S' = TyTf^Ty. We 
apply this method to the -F ipy topological superconductor. 

Dehn twists of px -F ipy superconductors. — We consider a 
Px -F ipy superconductor with the following Hamiltonian: 

H = — t ^ (clcr' + h.C.) — p, ^ c\.Cr 

irr’) 'f' 

( 2 ) 

+ X {^r.r'clcl +h.c.) 

{rr') 


Hereinafter we set the hopping amplitude to f = 1 as the 
energy unit, and is an odd-parity pairing 

which follows the px -F ipy pairing symmetry: Ar, 7 -i-e,i = 
A, Ar,r-i-e„ = iA. When \p\ < 4, the Hamiltonian describes 
a topological superconductor with Chern number ix = 1, also 
known** as the “weak pairing” phase; while in the “strong 
pairing” phase \p\ > 4 it describes a trivial superconductor. 
In the momentum space, the Bogoliubov-de Gennes (BdG) 
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Hamiltonian reads 


Hk =[— 2 (cos kx + cos ky) — fi]az 
+ A(sin kxCTx + sin ky<jy). 


(3) 


We notice that the px + ipy superconductor itself is a short- 
range entangled state, i.e. the topological entanglement en¬ 
tropy is zero. When put on a torus with periodic boundary 
conditions, there is a unique ground state and the S', T matri¬ 
ces are both scalars. 

However, since this is a fermionic system, one has four 
choices of boundary conditions for the fermions on a torus; 
Cr+L^ = e^'^'°Cr,Cr+Ly = Here ^x,y is the amount 

of magnetic flux threading through the holes of the torus 
and the distinct values are 0 and tt mod 27r due to the flux 
quantization. We notice that in a topological superfluid with 
gapless Goldstone modes, the phase 0 of the pairing or¬ 
der parameter(i.e. A = |A|e®®) has winding numbers Wx 
and Wy around the x and y directions of the torus due to 
the magnetic fluxes. Therefore, the ground state energy re¬ 
ceives the following contribution from the phase stiffness; 
Y / cPr (V0)^ = 27r^ps(Wx+Wy), where ps is the superfluid 
density per unit area. As a result, the four sectors correspond¬ 
ing to (wx,Wy) = ( 0 , 0 ), ( 1 , 0 ), ( 0 , 1 ), ( 1 , 1 ) are not degener¬ 
ate, and have the energy splitting independent of the system 
size-"^. To avoid this issue, one can couple the fermions to 
a dynamical 2D gauge field, then the Anderson-Higgs mech¬ 
anism will completely gap out the Goldstone mode. In our 
calculation, the order parameter is treated as a static extrinsic 
parameter, so the phase fluctuations do not enter. 

For Px + ipy superconductor on an Lx x Ly lattice, the 
four sectors split into two subspaces; three with even fermion 
parity and one with odd fermion parity. Exactly which one 
corresponds to the odd fermion parity depends on the sign of 
p and the parities of Lx and Ly. Interestingly, the p > 0 and 
p < 0 cases have distinct patterns of ground state fermion 
parities. If /r < 0, the ground state fermion parity does not 
depend on the parities of Lx and Ly, and the $3, = 0, <i>y = 0 
sector always has the odd fermion parity. 

We now consider implementation of a Dehn twist, by shear¬ 
ing the lattice along the y direction as depicted in Fig. 1(a). 
We weaken the bonds connecting (0, y) and (1, y), and at the 
same time strengthen the diagonal bonds connecting ( 0 , y) 
and (l,y -f 1). This is done by modifying the Hamiltonian 
H in Eq. (2) to H{X) = H -\- \Hi where 


^ (4) 

+ ~ +h.C.). 

V 

X going from 0 to 1 is one step of the Dehn twist. Repeat¬ 
ing this procedure for Ly steps accomplishes one full Dehn 
twist. The flux configurations are transformed from ($ 3 ,, ^y) 

to ($^-|-$y,$y). 


To have a well-defined Berry phase, we need to make sure 
that during the whole shearing process the system remains 
gapped (in the thermodynamic limit). To this end, notice 
that the modification Hamiltonian Hi basically creates a one¬ 
dimensional channel along the j/-direction, which can poten¬ 
tially localizes subgap states. More concretely, because during 
the entire twist the translation invariance along the y-direction 
is preserved, we can partially diagonalize the Hamiltonian by 
going to the Fourier space Cx{ky) = Ly^^"^ J2y 
and obtain a collection of one-dimensional Hamiltonians 
along X labeled by the momentum ky. The fluxe $j, merely 
change the quantization condition of the momentum. It is easy 
to see that for the purpose of studying subgap states, we only 
need to consider ky = tt and the ID Hamiltonian reads 

Hky=v = ~ ^ + h.c.) + {2 — p) 'y ( cl,Cx 

x^O X 

+ A^(c), 4 + 1 -f h.c.) (5) 

+ (2A - l)(4ci - Acgcj -f h.c.) 


Here Cx = Cx{tt). One immediately recognizes that this de¬ 
scribes a one-dimensional spinless p-wave superconductor, 
belonging to symmetry class D in the topological supercon¬ 
ductor periodic table.If 0 < /i < 4, the one-dimensional 
model is in the nontrivial topological phase. As the bound¬ 
ary coupling 2A — 1 changes sign, the ground state fermion 
parity of the one-dimensional Hamiltonian also changes and 
there must be a level crossing.As a result, gapless edge 
modes emerge along the twist line in the course of Dehn 
twist, which makes the evolution non-adiabatic in the ther¬ 
modynamic limit.On the other hand, when p < 0 Eq. (5) 
describes a trivial one-dimensional superconductor and there 
are no subgap modes for any value of A, so the bulk gap re¬ 
mains finite. For ky = 0, the boundary coupling (between the 
a; = 0 and a; = 1 sites) remains a constant, so there can not be 
any subgap states. In conclusion, the lattice Dehn twist is an 
adiabatic evolution only when p < 0. So in the following, we 
will set the chemical potential at p = —2, and study the Dehn 
twist on the px + ipy topological superconductor. 

S and T Matrices. — We now present numerical results of 
the T matrix for a px + ipy superconductor on a finite torus. 
We label the ground state in the flux sector ($ 3 ,, ) = 7 r(a, b) 

as with a,b = 0, 1. In this basis, the Ty matrix takes the 
following form; 


Ty = 


/ gi0oo 
0 
0 
0 


0 

0 

giOi 

0 


c.i©oi 

0 

0 


0 

0 

0 


( 6 ) 


Schematically, Ty\t^ab) = e‘^“*’| 4 'a+bmod 2 .&) transforms the 
fermion many-body state |'kab) to the new state | 4 'Q+bmod 2 ,b) 
evolved from the original state under a full Dehn twist along 
the y direction. 

In Fig. 2(a) we show the accumulation of Berry phase dur¬ 
ing the evolution. In each step of the twist, A is smoothly 
tuned from 0 to 1. For an Lx x Ly sized lattice, we need 
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FIG. 2; (a) The wave function overlap (upper panel) with the ref¬ 
erence state and the Berry connection d0/dA (lower panel) along 
the path of a full Dehn twist, with the initial flux configuration being 
(0,0) (in red) or (0, tt) (in green), (b) The Dehn twist Berry phase 
0 v.s. L^, fitted by 0 = a + j3L^ for large-L data to determine P 
(which is the same for all flux configurations). 
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FIG. 3: (a) The phase a v.s. system size L. q = 0 — piP is ob¬ 
tained from the calculated Berry phase 0 subtracted by pL^, where 
P determined from the fitting of the large-L data, (b) The jump of the 
universal phase a across the topological-trivial transition, fj, > —4 
(/r < —4) corresponds to the topological (trivial) phase. The cal¬ 
culation is performed on a 20 x 20 lattice. The sharp transition at 
^ = —4 is smoothed out by the finite-size effect. 


Ly steps of twists to complete a full Dehn twist. For exam¬ 
ple, the cases shown in Fig. 2(a) are of Ly — 9. Let |'I'(A)) 
be the many-body ground state of the Hamiltonian H{X). 
The Berry connection d0/dA can be calculated by compar¬ 
ing the phases of the nearby states |'F(A)) and |T'(A + dA)) as 
d0 = arg(T'(A)|'I'(A + dA)). To fix the gauge, we need to 
choose a reference state | thref) and make sure that there is no 
relative phase between |d>(A)) and I'll ref) along the way of the 
twist, i.e. arg(T'(A)|'kref) = 0 (see Appendix A for the calcu¬ 
lation of the overlap of two BCS wavefunctions and Appendix 
B for a more efficient algorithm in the momentum space). The 
upper panel of Fig. 2(a) shows the overlap of the ground states 
|T'(A)) with the reference state [dlref), which never vanishes, 
verifying that the Dehn twist is indeed performed adiabati- 
cally. The lower panel of Fig. 2(a) shows the Berry connec¬ 
tion through out the whole process, and the Berry phase can 
be obtained from the accumulation O — j> d0. 

It is worth pointing out that the Ty matrix is not completely 
diagonal; it leaves both the (0,0) and (7r,0) flux configu¬ 
rations invariant, but takes the (0, tt) configuration to (tt, tt) 
and vice versa. In order to get a well-defined Berry phase in 
the latter case, we complete the cycle by first starting from 


the (tt, tt) configuration and adiabatically turn off $ 3 , on the 
boundary bounds without twisting the lattice until we reach 
the ( 0 , tt) configuration, which is dubbed as a “flux removal” 
procedure; and then Dehn twisting the lattice back to the 
(tt, tt) configuration. Fig. 2(a) shows the whole evolution for 
both the (0, 0) case (which does not involves the flux removal) 
and the (0, tt) case (which involves the flux removal). It is evi¬ 
dent from the numerical result that the flux removal procedure 
does not give any contribution to the Berry phase. 

We carry out the calculation on a L x L sized lattice, and 
then perform a finite-size scaling of 0 to extract the universal 
phase a. We notice that for all the configurations, the values of 
P are almost identical (see Fig. 2(b)), so we can first determine 
P from the fitting of the large-L data, and then calculate a by 
subtracting PLP from the calculated 0, i.e. a = 0 — PLp. 
Fig. 3(a) shows that a approaches to its universal value in the 
thermodynamic limit as L — 00 . We notice that the value 
of P is actually related to the Hall viscosity of the fluid once 
we take the continuum limit: P = —rjH where rjH oi n for 
Px + ^Py superfluid with n being the number density**"’**’’ ’. 
We confirm numerically the linear dependence of /3 on n in 
the low-density regime (see Appendix C for more details). 

Up to numerical error, we And that aoo = aio = 7r/12 and 
aoi = ail = —7r/24. Plugging into Eq. (6), the T matrix is 
therefore 


Ty = e 


/eT 0 0 0 \ 

0 0 10 
0 10 0 
Vo 0 0 eT/ 


(7) 


Notice that we drop the overall non-universal phase 
The overall phase e 24 = e 24 ■ 2 is consistent with the chi¬ 
ral central charge c_ = 1/2. We can also perform the Dehn 
twist by shearing along the x direction, and calculate the cor¬ 
responding T matrix using the same method in the same set 
of basis. It is found that 


/e s' 

0 

0 

V 0 


0 

e~^ 

0 

0 


0 o\ 
0 0 
0 1 
1 oj 


( 8 ) 


Then the S matrix can be obtained from S = TyT~^Ty, be¬ 
cause the S transform corresponds to a 7r/2 rotation, which 
can be composed from three shear transformations (Dehn 
twists): 


,5 = TyT-% 


/eT 0 0 0\ 
0 0 0 1 
0 0 10 
V 0 100/ 


(9) 


If our model is viewed as fermions coupled to static Z 2 
gauge fields, the physical states must be gauge-invariant. So it 
is necessary to project to even fermion parity subspace which 
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corresponds to the lower 3x3 block. In the T^-diagonal basis, 


A 

0 

0 \ / 

0 

eT 

II 

0 

VO 

0 

-1/ \ 



Ty — e 


As expected, they match the topological T, S matrices of the 
Ising anyon model"’. It is interesting to observe that although 
the (0, 0) sector is thrown away in the gauge projection (be¬ 
cause its fermion parity is odd), the Berry phase of the Dehn 
twist is nevertheless still well-defined, which turns out to be 
the topological twist of the tt-Aux in the gauged model. 

TQFT Interpretation .— The above numerical results can be 
interpreted from a TQFT perspective. As we have mentioned, 
our model can be considered as a -f \py superconductor 
coupled to a static Z 2 gauge field, and the TQFT description 
is an Ising anyon model , which has three topological charges, 
1, (T and 'ijj. Here ip corresponds to the fermion in the lattice 
model, and a corresponds to the tt flux with a localized Majo- 
rana zero mode. They satisfy the well-known fusion rules: 


axa = l + ip,axip = a. 


( 11 ) 


On a torus, the topologically-protected degenerate ground 
states can be labeled by the topological charge that going 
through one of the two holes, e.g. the longitudinal one (often 
refered to as an anyon flux), and a topological charge mea¬ 
surement can be performed along the meridian cycle to distin¬ 
guish the different ground states. The S and T matrices can 
be understood as basis transformations on the torus. 

For the + ipy superconductor, the three ground states 
with even fermion parity are unitarily related to the states 
with definite topological charges l^a^ip. However, the odd 
fermion parity sector (0,0) requires a different interpretation: 
formally this state does not exist on a closed surface since it 
carries an odd Z 2 charge. Instead one has to use a punctured 
torus with a ip charge line coming out of the punctured hole, 
see Fig. 4(a). Judging from the fusion rule, this means that 
there must be a cr charge line going in the longitudinal direc¬ 
tion (see Fig. 4(a) for an illustration), which explains why the 
Berry phase of the Dehn twist on the (0,0) sector is equal 
to the topological twist of a. The S matrix is now given by 
as shown in Fig. 4(b), and our numerical result indeed 
agrees with the diagrammatic result. We can also understand 
this result intuitively as follows: assuming there is SO(2) ro¬ 
tation symmetry, a single fermion acquires — 1 phase under a 
full 27r rotation(i.e. the topological spin). Therefore a | rota¬ 
tion would results in phase, which is consistent with the 
numerical calculation. 

Topological Phase Transition .— The Hamiltonian Eq. (2) 
can also describe a trivial superconductor if the chemical po¬ 
tential p lies outside the band, i.e. \p\ > 4, separated from the 
topological phase by a quantum phase transition. In the above 
calculation, we have always set p = —2 in the topological 
phase. Now we tune the chemical potential across the phase 
transition, and perform the same calculation of the Dehn twist 
Berry phase. Fig.3(b) shows that the universal phases a iden¬ 
tically drop to zero across the topological-trivial transition as 





FIG. 4: (a) The odd fermion parity sector can be view as a fermion 
charge line ip going out from the punctured torus, (b) Diagrammatic 
definition of the S matrix in the odd fermion parity sector. 


the chemical potential is tuned below —4. So if we start with 
the trivial superconductor, the S and T matrices will be 


Ty = 


a 0 0 ON 

0010 
0100 
vO 0 0 1 / 




/I 0 0 ON 
0 0 0 1 
0 0 10 
\o 1 0 oy 


( 12 ) 


For the trivial superconductor, all the flux configurations lead 
to even fermion parity states, so the S and T matrices remain 
the same as Eq. (12) after gauging the fermion parity, which 
can be transformedto the standard form of the S and T 
matrices of the Z 2 toric code model. The TQET description 
of a Z 2 gauged trivial superconductor is indeed a toric code 
theory, which has four topological charges 1, e, m and ip, 
whose topological twists are given by the eigenvalues of Ty-. 
1 ,1,1, —1. Ip corresponds to the fermion in the lattice model, 
m corresponds to the tt flux in the trivial superconductor, and 
e is the bound state of ip and m. 

TSC with Higher Chern Number .— A topological super¬ 
conductor with Chern number v > 0 can be considered as 11 
copies of px + ipy superconductors. Each copy will accumu¬ 
late the same amount of Berry phase during the Dehn twist, so 
the S and T matrices simply follow from Eq. (7) and Eq. (9) 
with all the Berry phases multiplied by 1 / times, Eor odd v, 
the even fermion parity sector corresponds to the lower 3x3 
block; while for even v, all the four basis states are of even 
fermion parity. 

Conclusions and Discussions .— In this work, we calcu¬ 
lated both the S and T modular matrices of two-dimensional 
topologically ordered phases that can be constructed by gaug¬ 
ing the fermion parity from a free-fermion topological super¬ 
conductors. The modular matrices can be obtained as the non- 
Abelian Berry phases among the ground states in different 
gauge flux sectors. We use Dehn twists along both the merid¬ 
ian and the longitudinal cycles of the torus as the generators of 
the modular group, which can be conveniently implemented 
by smoothly shearing the lattice in both directions. It will be 
interesting to generalize this method to calculate the modu¬ 
lar data of interacting topological phases on lattices such as 
fractional Chern insulatorsas well as G-crossed modular 
data for G-symmetry-protected phases^*, which will be left 
for future investigations. 
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Appendix A: Overlap of BCS wavefunctions 


Appendix B: Dehn Twist in Momentum Space 


In this appendix we summarize the algorithm to calcu¬ 
late the overlap of BCS wavefunctions. Let us rewrite the 
complex fermions Cr in terms of Majorana fermions Xz, s.t. 
Cr = Xr,i + iXr', 2 - In general the BCS Hamiltonian can be 
written in the Majorana basis 

ij 

where A is an external parameter. Diagonalizing the Hermitian 
matrix iA, one finds its eigen states along with the eigen 
energies i^lun = EnUn- We only consider the gapped 
fermion system, so the spectrum contains half of the states 
Un with negative energy < 0. The fermionic many-body 
ground state jT*) of H is simply constructed by occupying all 
the negative energy single-particle states. Let us arrange the 
negative energy states Un (in the form of column vectors) into 
a matrix 


(A2) 

such that each column of (7 is a state vector Un- Of course, this 
matrix (7(A) also depends on the external parameter A. Then 
the overlap of two BCS ground states |T'(Ai)) and |T'(A 2 )) 
can be expressed in terms of the corresponding U matrices as 


The lattice Hamiltonian in Eq. (2) can be easily diagonal¬ 
ized in the momentum space as Eq. (3). On the twisted lattice, 
the translational symmetry that is perpendicular to the twist 
should be understood as the translation followed by a shearing 
along the Dehn twist direction, therefore the quasi-momenta 
(eigenvalues of the translation symmetry) are quantized to (on 
an L X L square lattice) 


^ f 1 0\ / 2TTnx + 

L \-Ty V + ^yj ’ 


'^xi^y ^ ^7 (HI) 


where Ty evolves from 0 to 1 corresponding to a full Dehn 
twist along the y-direction, and denotes the initial 

flux configuration (at the beginning Tj, = 0 of the Dehn twist). 

Eor each k in the Brillouin zone, we can simply diagonalize 
the 2x2 single-particle Hamiltonian given in Eq. (3) to 
obtain the Bloch state vector Uk, i.e. HkUk — E^Uk (with 
Ek < 0). Then each occupied single-particle state is labeled 
by both its quasi-momentum k and the corresponding Bloch 
state vector Uk, denoted as |fc, Uk)- The overlap between two 
single-particle states is given by 

{k,Uk\k ,Uk') = r n,, -(B2) 

Lsin — k')\ 


(T'(Ai)|T-(A 2 )) = Pf [C/I(Ai)C/(A 2 )] . (A3) 

With this, one can calculate the overlap between any pair of 
BSC ground states. 

However the eigen states obtained from diagonalization 
of A is subject to an overall phase degree of freedom, and 
hence the U matrix is also ambiguous. To remove this phase 
ambiguity and fixed the gauge, we need to choose a reference 
state I T'ref), so that we can define the Berry phase accumulated 
from A to A -f dA to be 

d0=arg(vI/(A)|T'(A + dA)) 

- arg(T'(A)|T'ref) + arg(T'(A -f dA)|T'ref)- 

Alternatively, one may first rectify the overall phase of |'I'(A)) 
with respect to the reference state, such that arg('I'(A) jT'ref) = 
Oforall A, then Eq. (A4) is reduced to d0 = arg(T'(A)|T'(A-|- 
dA)). On every given small segment dA, d0 may vary with 
the choice of the reference state | ■ But after integrating 

over a closed path, the Berry phase & = j> d0 is well-defined 
and independent of I'l'ref)- In the calculation of the Berry 
phase, it is import to check that the overlap to the reference 
state |(T'(A)|'kref)| never becomes vanishingly small, other¬ 
wise the result will suffer from severe numerical error. 

In our calculation, for ctoo and aio, the reference state is 
chosen to be the initial state Ithref) = |dt(A = 0)). While for 
aoi and an, the reference state is chosen to be the ground 
state of H with the boundary condition being open in the x- 
direction and anti-periodic in the y-direction. 


With this, the overlap of two fermionic many-body states that 
are parametrized by Ty i and 2 follows from Eq. (A3) as 

(^('ry.l)l^(T'y, 2 )) = Pf[r„,i(fe,Mfc|fc',Ufe')T„, 2 ], (B3) 

where \k,Uk)Ty denotes the single-particle state that is ob¬ 
tained from the Hamiltonian at the parameter Ty, and the Pfaf- 
flan matrix index runs over k (row index) and k' (column in¬ 
dex) through out the Brillouin zone. Then the calculation of 
the Berry phase is similar to the real space algorithm given 
in Eq. (A4). In this momentum space representation, the flux 
insertion/removal procedure can be simply implemented by 
smoothly tuning (4)a;, 4>y). The computation efficiency can be 
improved significantly using the momentum space formalism, 
because the single-particle Hamiltonian is very easy to diago¬ 
nalize in the momentum space. 


Appendix C: Hall viscosity from the Lattice Model 

We change the number density of the + ipy superfluid 
by tuning the chemical potential p in the model Hamiltonian 
Eq. (2). Eor each density fi, we calculated the Dehn twist 
Berry phase 0 = a -I- + 0{L~‘^) for various L, and 

extract the coefficient /3. /3 is related to the Hall viscosity rju 
of the fluid in the continuum limit: f3 = —rjH which is pro¬ 
portional to the number density n. It is verified that (3 indeed 
scales linearly with h in the dilute limit n —0, as expected. 



6 
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FIG. 5: Hall viscosity scales with the number density n linearly in 
the dilute limit n —>■ 0. 
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Here we notice a potentially confusing convention in defining the 
boundary conditions of fermions. If we actually consider fermions 
living on a toms which embeds in R®, the boundary conditions of 
Cr is shifted relatively to the actual flux by rr, because of the spin 
connections on the toms. 

This is the essence of the drr Josephson effect of class D topolog¬ 
ical superconductors. 

Even in a finite system where one can suitably chose Ly to avoid 
the ky = TT point and maintain a finite-size gap during the whole 
twist, the subgap mode still has significant contribution to the 
Berry phase. 




